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ABSTRACT 

RAPID  INITIALIZATION  OF  INERTIAL  NAVIGATION  SYSTEMS 
THROUGH  PARAMETER  ESTIMATION 

The  accuracy  of  an  aircraft  inertial  navigation  system  de» 
pends  upon  the  accuracy  with  which  the  system  is  initially 
aligned.  One  procedure  for  initial  alignment  involves  the  uaa 
of  an  external  reference.  This  method  utilizes  equipment  which 
is  much  too  elaborate  for  normal  operational  use.  An  alternate 
procedure  uses  the  system's  inertial  sensors  in  a  self-contained 
method.  If  sufficient  time  were  available,  the  self-contained 
method  could  achieve  accuracies  commensurate  with  the  sensor 
accuracies}  however,  in  an  operational  environment  it  is  usually 
necessary  to  sacrifice  some  accuracy  in  the  interest  of  achieving 
a  more  rapid  initiation.  This  dissertation  investigates  the 
methods  presently  available  for  initialization  of  an  inertial 
platform  in  an  azimuth  wander  or  free  azimuth  instrumentation^  sad 
presents  a  new  method  for  rapid  initialization.  The  paramount 
problem  is  the  determination  of  the  initial  azimuth  angle  in 
minimum  time  in  the  presence  of  random  gyro  drifts,  random  aoeol- 
erometer  drifts,  and  measurement  ncise. 

A  linear  error  model  of  the  inertial  platform  is  developed. 
The  model  contains  all  significant  cross-coupling  terms.  The 
inertial  component  random  drifts  are  modeled  as  time  correlated 
random  processes  and  the  measurement  noise  is  represented  as  a 
white  Gaussian  process.  State  space  equations  for  the  error 
model  are  then  formulated.  The  problem  of  determining  the 
initial  azimuth  wander  angle  is  then  identified  as  a  parameter 


eatimation  problem  where  the  parameter  can  assume  any  of  a  con¬ 
tinuum  of  values  (from  0  to  2iy).  The  two  methods  of  solving  pa¬ 
rameter  estimation  problems  currently  available  in  the  literature 
ere  presented.  The  first  method  allows  parameter  estimation 
when  the  parameter  can  assume  a  continuum  of  values;  however,  the 
method  is  not  time  optimal.  The  second  method  examined  is  time 
optimal;  however,  it  is  constrained  to  problems  where  the  param¬ 
eter  can  assume  only  a  finite  number  of  possible  values.  The 
second  method  is  then  extended  to  permit  time  optimal  parameter 
estimation  where  the  parameter  can  assume  a  continuum  of  values. 

The  parameter  estimation  method  developed  utilizes  an  array 
of  minimum  variance  filters.  Each  filter  (referred  to  as  an 
elemental  filter)  is  initialized  with  an  estimate  of  the  unknown 
parameter.  One  element  of  the  filter  state  vector  is  related  to 
the  parameter  such  thet  feedback  can  be  used  to  continually  up¬ 
date  the  estimate  of  the  parameter.  The  elemental  filters'  pa¬ 
rameter  values  are  then  combined  to  form  the  overall  parameter 
estimate.  A  simple  weighting  scheme  is  used  in  the  combining 
procedure.  A  variance  term  for  the  parameter  estimate  is  also 
computed  so  that  the  initialization  procedure  can  be  terminated 
when  a  predetermined  variance  is  achieved. 

The  equations  for  a  discrete  minimum  variance  filter  with 
internal  feedback  ere  then  presented.  A  computer  algorithm  is 
developed  for  the  elemental  filter  of  the  parameter  estimator. 

A  method  of  applying  the  parameter  estimation  technique  to 
determine  the  initial  azimuth  wander  angle  is  then  formulated. 


iii 


The  platform  controller  in  then  developed  end  the  overall  system 
described , 

A  computer  simulation  of  the  rapid  initialization  of  an 
aeimuth  wander  system  through  the  use  of  the  parameter  estimation 
technique  is  discussed.  Results  of  the  simulation  are  presented 
showing  reduction  of  an  initial  wander  angle  error  variance  from 
(1.5  degrees)3  to  (6  minutes  of  arc)3  after  approximately  threa 
minutes  of  real  tiifia.  This  represents  a  twofold  improvement  in 
initialization  time  over  a  state-of-the-art  system  presently 
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I.  Introduction 

The  accuracy  of  an  aircraft  inertial  navigation  system  de¬ 
pends  upon  the  accuracy  with  which  the  system  is  initially  aligned,1 
One  procedure  for  initial  alignment  involves  the  use  of  an  ex¬ 
ternal  reference.  This  method  utilizes  equipment  which  is  much 
too  elaborate  for  normal  operational  use.  An  alternate  procedure 
uses  the  system's  inertial  sensors  in  a  self  contained  method. 

If  sufficient  time  were  available,  the  self  contained  method 
could  achieve  accuracies  commensurate  with  the  sensor  accuracies} 
however,  in  an  operational  environment  it  is  usually  necessary  to 
sacrifice  some  accuracy  in  the  interest  of  achieving  a  more  rapid 
initiation.  In  this  dissertation  the  problem  of  determining  the 
initial  orientation  of  an  inertial  platform  is  investigated. 

An  inertial  navigation  system  consists  of  three  sub-systems: 
the  inertial  measurement  unit  or  sensor  package,  the  computer  in 
which  the  navigation  equations  are  solved,  and  the  display. 

Prior  to  operation  as  an  autonavigator,  two  sets  of  quantities 
are  required  as  initial  conditions  fcr  the  inertial  navigation 
system  equations. 

First,  the  true  initial  position  of  the  inertial  measurement 
unit  in  terms  of  a  navigation  system  coordinate  frame  must  be 
accurately  known.  In  addition,  the  angular  orientation  of  the 
axes  of  the  inertial  measurement  unit  with  respect  to  the  naviga¬ 
tion  system  coordinate  frame  must  also  be  accurately  known.  In 
this  investigation  it  is  assumed  that  initial  position  is  known  a 
priori  and  the  angular  orientation  is  to  be  determined. 
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At  this  point  it  is  well  to  define  the  various  coordinate 
frames  which  will  be  referred  to  in  the  ensuing  paragraphs. 

First,  an  inertial  frame  with  its  origin  at  the  center  of  die 
earth  will  be  the  basic  fixed  frame.  The  navigation  reference 
frame  will  be  the  conventional  latitude-longitude  frame,  which  is 
also  earth  centered  but  is  rotating  with  respect  to  the  fixed 
frame.  At  a  point  in  this  navigation  frame  (a  specific  latitude- 
longitude  point  on  the  earth's  surface)  a  locally  level  triad  can 
be  defined  by  Xjj,  Yn,  Zn  where  Xn  is  level  and  pointing  east,  Yn 
is  level  and  pointing  north,  and  2^  is  up.  Finally,  the  inertial 
measurement  unit  will  have  imbedded  a  set  of  axes  Xp,  Y^,  Z& 
which  are  mutually  orthogonal.  The  angular  orientation  of  these 
Inertial  measurement  unit  axes  with  respect  to  the  local  level 
triad  must  be  accurately  known  for  proper  initialization. 

The  inertial  measurement  unit  in  this  investigation  is  the 
usual  Schuler- tuned  gimballed  platform  for  terrestrial  navigation.2 
Within  the  platform  assembly  are  mounted  thrae  gyros,  the  sensi¬ 
tive  axes  of  which  form  a  mutually  orthogonal  triad,  and  two 
accelerometers.  When  in  operation,  two  of  the  sensitive  gyro  axes 
are  oriented  in  the  horizontal  or  locally  level  plane  and  the 
sensitive  axes  of  the  two  accelerometers  are  collinear  with  these 
gyros'  sensitive  axes.  The  third  gyro  has  its  sensitive  axis  in 
the  vertical  direction  and  is  referred  to  as  the  azimuth  gyro. 

With  respect  to  the  local  level  reference  triad,  the  usual  con¬ 
vention  is  to  have  one  of  the  level  gyros'  sensitive  axes  pointed 
east  and  the  other  pointed  north.  These  gyros  are  oometimes 


referred  to  ss  the  east  gyro  and  north  gyro  respectively,  and  de¬ 
fine  the  Xp  and  Yp  axes  of  the  inertial  measurement  unit.  “Hie 
azimuth  gyro  sensitive  axis  (pointed  up)  defines  the  inertial 
measurement  unit  Zp  axis. 

For  a  stationary  system  (that  is,  no  vehicla  velocity),  one 
method  for  obtaining  the  initial  angular  orientation  of  the  level 
Xp,  Yp  platform  axes  with  respect  to  the  Xn,  Yn  axes  is  through 
use  of  an  external  true  north  reference.  The  orientation  of  the 
reference  in  the  local  level  axes  must  be  fixed  and  accurately 
known.  The  information  transfer  from  the  external  reference  is 
accomplished  through  use  of  precision  electronic  theodolites  and 
mirrors  located  on  the  platform.  This  method  of  initiation  is 
limited  to  operational  environments  which  can  support  this  cm- 
bersome  alignment  procedure.  For  most  land  based  tactical  air¬ 
craft  a  self  contained  alignment  capability  is  necessary  to  pro¬ 
vide  operational  flexibility. 

The  primary  method  for  self  contained  alignment  of  an  inertial 
platform  is  commonly  referred  to  as  gyrocompassing.®  For  a 
stationary  system,  gyrocorapassing  consists  of  two  phases:  platform 
leveling  and  azimuth  alignment.  The  simplest  form  of  leveling 
consists  of  feeding  the  accelerometer  gravity  measurements  to  the 
appropriate  gyro  torquers.  The  platform  is  then  torqued  until 
the  sensitive  axes  of  the  level  gyros  are  in  the  locally  level 
plane.  The  control  loop  gains  and  compensation  networks  are  set 
for  the  desired  response.4  The  more  complex  azimuth  alignment 


phase  is  described  in  the  following  paragraphs. 

To  maintain  the  platform  orientation  in  the  local  level  frame 
it  is  necessary  to  tcrque  the  platform  axes  at  the  same  rate  that 
the  local  level  frame  is  moving  with  respect  to  the  basic  fixed 
frame.  Otherwise,  the  platform  would  maintain  its  orientation  in 
the  inertial  or  fixed  frame.  If  the  platform  Xp  and  axes  are 
oriented  east  and  north  in  the  local  level  frame,  a  torquing  rate 
must  be  applied  to  the  north  pointing  gyro  to  maintain  the  plat¬ 
form  in  a  level  position.  No  torquing  of  the  east  pointing  gyro 
is  necessary  since  its  sensitive  axis  io  orthogonal  to  the  rota¬ 
tion. 

If  the  Xp  and  Yp  platform  axes  are  not  aligned  to  east  (Xjj) 
and  north  (Yn)  in  the  local  level  frame,  the  platform  will  appear 
to  be  tilting  when  viewed  in  this  frame.  This  tilting  property 
is  used  to  align  the  Xp  and  Yp  platform  axes  with  the  X^  and 
local  level  axes  in  the  azimuth  alignment  phase  of  gyro compassing* 

When  the  platform  Xp  axis  is  not  coincident  with  in  the 
local  level  frame,  a  rotation  about  the  Xp  axis  as  viewed  in  the 
local  level  frame  will  cause  the  north  accelerometer  to  sense  a 
component  of  gravity.  This  signal  is  then  applied  to  the  Zp  or 
azimuth  gyro  torquer.  The  platform  is  then  rotated  until  the 
Xp  axis  is  brought  into  coincidence  with  the  local  level  X„  axis 
and  gyrocompassing  completed. 

In  order  to  maintain  the  platform  level,  it  is  necessary  to 
apply  torquing  signals  to  both  the  north  gyro  and  the  azimuth  gyro 
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at  latitudes  other  than  the  equator.  The  vertical  component 
appears  because ,  when  viewed  from  a  point  on  the  earth ' s  surface , 
the  earth^s. rotation  vector  can  be  resolved  into  a  horizontal  or 
level  component  and  a  vertical  component.  This  level  component 
•is  greatest  at  the  equator  and  zero  at  the  poles,  varying  in  a 
sinusoidal  manner.  The  azimuth  alignment  phase  of  gyrecoinpassing 
is  more  difficult  in  the  polar  regions  due  to  the  smaller  level 
component  of  the  earth's  rotation.8 

Some  inertial  platforms  in  use  today  are  not  aligned  in 
azimuth.  Instead,  the  angle  between  the  platform  Yp  axis  and 
local  level  Yn  is  determined.  This  angle  is  generally  referred 
to  as  the  "azimuth  wander"  angle  and  the  inertial  system  implemen¬ 
tation  as  an  azimuth  wander  mechanization .*» 7 

In  an  azimuth  wander  mechanization,  torquing  signals,  which 
maintain  the  platform  level,  are  applied  to  both  the  Yp  and  Xp 
gyros.  The  torquing  signals  are  determined  by  using  the  azimuth 
wandar  angle  to  revive  the  level  angulor  earth's  rate  into 
components  along  the  actual  orientation  of  the  level  platform  axes. 
The  vertical  component  of  the  earth's  angular  rate  is  integrated 
and  the  azimuth  wander  angle  continuously  corrected.  This 
dissertation  will  consider  a  method  for  rapid  determination  of 
the  initial  azimuth  wander  angle. 

Initial  alignment  of  a  locally  level  azimuth  wander  inertial 
navigator  is  normally  accomplished  in  either  of  the  two  ways 
previously  discussed,  that  is,  through  use  of  an  external  reference 
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or  a  modified  fora  of  gyrocompassing.  When  an  external  reference 
in  utilized  the  platfora  is  first  leveled ,  and  then  the  actual 
initial,  azimuth  wander  angle  determined  through  uae  of  optical 
instruments .  In  modified  gyrocompassing ,  a  coarse  estimate  of 
the  azimuth  wander  angle  is  initially  obtained  from  an  independent 
reference,  such  as  a  magnetic  compass.  The  platform  is  then 
leveled  and  this  estimate  of  the  azimuth  wander  angle  used  for  the 
computation  of  the  angular  rates  to  be  applied  to  the  level  gyro 
torquers.  If  tilting  (with  respect  to  the  local  level  exes) 
occurs,  the  wander  angle  is  updated.  The  updated  value  is  based 
on  accelerometer  gravity  measurements .  The  platform  is  then  re- 
leveled  and  new  torquing  rates  computed  for  the  level  gyro  tor¬ 
quers.  This  sequence  is  repeated  untii  the  accelerometer  measure¬ 
ments  remain  below  a  pro-determined  level  for  a  fixed  attount  of 
ti-xe . 

Self  contained  alignment  schemes  generally  require  between  • 
twenty  and  thirty  minutes  of  time  to  achieve  an  initial  azimuth 
accuracy  on  the  order  of  six  minutes  of  arc.8  Three  factors 
which  affect  alignment  time  and  accuracy  are  inertial  component 
warm-up  behavior,  random  inertial  component  drifts,  and  measure¬ 
ment  noise.9*10  During  the  past  few  years,  significant  accomplish¬ 
ments  in  the  area  of  thermal  modeling  of  inertial  component  drifts 
during  warm-up  have  been  reported11* 12 ,  and  the  availability  of 
the  digital  airborne  computer  has  made  thermal  compensation 
feasible.  This  dissertation  presents  a  method  of  initiating  a 


thermally  compensated  azimuth  wander  inertial  system  in  minimum 
time  (less  than  five  minutes). 

The  subsequent  paragraphs  present  the  organization  and  con¬ 
tents  of  this  dissertation. 

First,  a  linear  error  model  of  the  inertial  platform  in  the 
azimuth  wander  mechanization  is  developed.  The  model  contains 
all  significant  cross-coupling  terms.  The  inertial  component 
random  drifts  are  modeled  as  time  correlated  .random  processes  and 
the  measurement  noise  is  represented  as  a  white  Gaussian  process. 

Next,  state  space  equations  for  the  error  .model  are  formu¬ 
lated.  The  problem  of  determining  the  initial  azimuth  wander 
angle  is  then  identified  as  a  parameter  estimation  problem  where 
the  parameter  can  assume  any  of  a  continuum  of  values  (from  C  to 
2ir). 

The  two  methods  of  solving  parameter  estimation  problem  a 
currently  available  in  the  literature  are  presented.  The  first 
method  allows  parameter  estimation  when  the  parameter  can  assume 
a  continuum  of  values^  however,  the  method  is  not  time  opitmal. 
The-  second  method  examined  is  time  optimal,  however,  it  is  con¬ 
strained  to  problems  where  the  parameter  can  assume  only  a  finite 
number  of  possible  values.  The  second  method  is  then  ex tended 
to  permit  time  optimal  parameter  estimation  where  the  parameter 
can  assume  a  continuum  of  values. 

The  parameter  estimation  method  developed  utilizes  an  array 
of  minimum  variance  filters.  Each  filter  (referred  to  as  an 
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elemental  filter)  is  initialized  with  an  estimate  of  the  unknown 
parameter.  One  element  of  the  filter  state  vector  ie  related  -to 
the  parameter  such  that  feedback  can  be  used  to  continually  update 
the  estimate  of  the  parameter.  The  elemental  filters'  parameter 
values  are  then  combined  to  form  the  overall  parameter  estimate. 

A  simple  weighting  scheme  is  U3ed  in  the  combining  procedure.  A 
variance  term  for  the  parameter  estimate  is  also  computed  sc  that 
the  initialization  procedure  can  be  terminated  when  a  predetermined 
variance  is  achieved. 

The  equations  fop  a  discrete  minimtot  variance  filter  with 
internal  feedback  ar*  then  presented.  A  computer  algorithm  is 
developed  for  the  elemental  filter  of  the  parameter  estimator. 

A  method  of  applying  the  parameter  estimation  technique  to 
determine  the  initial  oziauth  wander  -angle  is  thsn  formulated. 

The  platform  controller  is  developed  and  the  overall  systwt 
described. 

A  computer  simulation  of  the  rapid  initialisation  cf  an 
azimuth  wander  system  through  the  use  of  the  parameter  estimation 
v9chnici%t  is  discussed.  Results  of  the  simulation  are  presented 
showing  reduction  of  an  initial  wander  angle  error  variance  from 
(l.S  degrees;2  to  1,6  minutes  of  arc)2  after  approximately  three 
Juicures  of  real  time. 

Conclusions  of  the  investigation  are  then  presented. 
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II.  The  System  Error  Model 

A«  stated  in  the  introduction,  the  distinguishing  character¬ 
istic  between  the  conventional  Schuler- tuned  gimballed  inertial 
system  and  an  azimuth  wander  system  is  the  absence  of  a  torquer 
on  the  vertical  or  azimuth  gyro.  Both  systems  normally  operate 
by'  keeping  their  platforms  locally  level  and  use  the  rotating 
navigation  coordinate  .frame  as  a  geographic  reference.  Due  tc 
the  absence  of  the  azimuth  gyro  torquar  in  the  azimuth  wander 
mechanization,  the  arbitrary  azimuth  angle  must  be  precisely 
known.  The  azimuth  wander  angle  is  used  in  a  rotation  matrix  to 
transform  platform  measurements  from  the  inertial  measurement 
unit  axes  {Xp,  Yp,  Zp)  to  the  local  level  navigation  axes  (X^ 

V  Vi- 

Due  to  the  earth's  rotation  the  navigation  reference  fxtme 
rotates  about  the  fixed  inertial  frame.  This  rate  cf  rotation  is 
referred  to  as  the  earth's  rate  and  will  be  represented  by  &  . 

For  a  point- on  the  navigation  axes,  the  level  and  vertical  com¬ 
ponents  of  il  can  be  computed  by 

ilyn  =  ill' COS  X  2-1 

o 

Qxr  —  Xl-ain  X  2-2 

where  X  is  the  latitude  of  the  point  and  the  subscripts  xn’ 
and  Z-j  represent  directions  along  the  local  level  triad.  Figure 
2-1  illustrates  the  resolution  of  the  earth's  rate  into  horizontal 
and  vertical  components  at  a  point,  ?,  on  the  earth's  surface. 
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The  platform  axes  (Xp3  Yp,  and  Zp)  tend  to  remain  fixed  with 
respect  to  the  fixed  inertial  frame.  To  maintain  the  platform 
axes  fixed  in  the  navigation  frame  it  is  necessary  to  apply  a 
signal  to  the  gyro  torquers  which  Vi 11  cause  the  platform  to 
process  at  the  angular  rates  computed  in  2-1  and  2-2.  If  the 
platform  is  initially  level  and  Yp  pointing  north,  the  application 
of  these  angular  torguing  rates  trill  maintain  the  platform  level 
and  north  pointing. 

Per  an  azimuth  wander  mechanization  the  orthogonal  platform 
level  axes  (Xp,  Yp)  will  be  oriented  arbitrarily  with  respect  to 
the  local  level  navigation  axes.  The  azimuth  wander  angle,  <tx, 
is  defined  as  the  angle  which  the  local  level  axes  must  he  rotated 
through:  in  a  counter-clockwise  direction,  to  become  coincident 
with  the  plavfom  level  exes.  Figure  2—2  illustrates  an  azimuth 
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0Z  (t) 


d  v 


+  $z  (0) 


2-6 


Hence,  the  value  of  02  in  the  rotation  matrix  must  be  continuously 
updated. 

The  rotation  matrix  in  equation  2-5  is  extremely  important 
in  an  azimuth  wander  mechanization.  Aside  from  resolving  level 
earth's  rate  into  components  along  the  platform  axes,  the  inverse 
of  this  matrix  is  the  transformation  from  the  level  platform 
axes  to  the  level  navigation  axes.  Therefore,  the  platform 
acceleration  measurements  end  platform  velocities  must  be  operated 
on  by  this  transformation  matrix  prior  to  use  in  the  navigation 
mode. 


After  platform  leveling,  an  initial  estimate  of  0z  is  used 

A 

in  the  rotation  matrix.  This  estimate  is  designated  0Z,  and  the 
difference  between  0Z  and  0Z  is 

A  0Z  =  $z  *>  2-7 

where  A  0Z  is  the  azimuth  wander  error  angle. 

Platform  level  axes  torquing  rates  are  then  computed  using 
equation  2-5  with  0Z  replaced  by  ^z.  The  computed  torquing 
signals  are  then  applied  to  the  level  gyro's  torquers.  The  error 
in  torquing  for  each  level  gyro  due  to  the  difference  between 


0Z  and  0Z  is 


=  iiyn«  cos  $z  -  Cl  yn»  cos  0Z 


E„  =  Hv.-sin  0,  -  y  » sin  '0„ 


2-8 

2-9 


y 

-x  =  ftyn,sin  '“z  -  “yn' 

where  Ey  is  the  error  in  earth  rate's  torquing  to  the  Y  gyro,  and 

Ex  is  the  error  in  torquing  to  the  X  gyro. 
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Using  equation  2-7  to 
and  2-9  results  in 

By  -  ft  yn*cos-  (0Z  - 
Bx  =  ft  yn*  sin  (#z  - 
Using  siaple  trigonometric  identities  for  the  expressions  with 
angle  difference  yields 

By  =  fit  yn  Qcos  cos  A0Z  +  sin  0Z  sin  A0Z  -  cos  0^2-12 
Ex  =  ftyn  C  sin  #z  cos  A  #z  -  cos  #2  sin  A0Z  -  sin  0^j2“13 
An  initial  coarse  estimate  of  0Z  can  be  obtained  from  a 
magnetic  compass.  After  correction  for  magnetic  variation,  the 
azimuth  wander  error  angle,  A#z,  is  of  the  order  of  1.5  degrees 
(standard  deviation).  At  the  conclusion  of  the  initiation  procoss 
an  error  of  0.1  degrees  (standard  deviation)  is  desirable  to 
minimiza  the  navigational  error  due  to  incorrect  initial  azimuth 
information. 1 

The  azimuth  wander  error  angle,  A  0_,  will  therefore  range 
in  value  from  an  initial  worst  case  value  of  approximately  4.5 
degrees  (three  standard  deviations)  to  approximately  0.1  degrees 
at  the  conclusion  of  the  initiation  process.  Under  these  cir¬ 
cumstances  the  following  small  angle  approximations  can  be  used 
cos  A  #z  =  * 

2-14 

sin  A  =  A  0Z 

Eq-iations  2-12  and  2-13  can  now  be  rewritten  as 

By  =  ft  yn  •  sin  #2  *  A  2-15 

Ex  =  “  ft  yn  •  cos  02  *  A  0Z 


substitute  for  0Z  in  equations  2-8 

A  0Z)  -  ftynaccs0z  2-10 

A  0Z )  -  ft  yn»  sin  0Z  2-11 


2-16 
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These  torquing  errors  cause  the  platform  to  tilt  in  the  navigation 

frame.  Let  A  #  and  A  0  represent  the  angles  between  the 
x  y 

locally  level  plane  and  the  platform's  Yp  end  Xp  axes.  That  is, 

A  0X  is  the  angle  between  Yp  and  the  locally  level  plane 
measured  by  a  right-handed  rotation  about  Xp,  and  A  0y  is  the 
angle  between  Xp  and  the  locally  level  plane  measured  by  a  right- 
handed  rotation  about  Yp.  After  the  initial  leveling,  these 
angles  can  be  considered  zero;  however,  the  application  of  the 
earth's  rate  torquing  signals  with  an  error  of  A  0Z  causes  a 


tilt  which  can  be  computed  from 


AeL(t) 


By  dt  = y  ft  yn  •  si 

/*  E*  at  *  -/*  a  yn  • 


sin  0Z*  A  0Z  dt 


2-17 


cos  GL,  *  A  0_dt  2-18 

25  2 


These  tilts  will  cause  the  level  axes'  accelerometers  to 
sense  gravity  components.  The  values  of  these  gravity  components 


will  be 


Ay  =  g*sin  A  0X  2-19 

Ax  -  -g#sin  A  2-20 

where.  Ay  and  Ax  are  the  platform  level  axes'  accelerometer  measure¬ 
ments.  The  minus  sign  in  2-20  results  from  assuming  a  right- 
handed  convention  for  the  platform  axes  and  a  gravity  vector  (g) 
in  the  positive  vertical  direction.  Therefore,  a  positive  error 
torque  to  the  platform  Xp  gyro  will  cause  the  Y  gyro  axis  to  lift. 

A  positive  g  component  will  then  be  sensed  by  the  Yp  accelerometer. 
Similarly,  a  positive  error  torque  to  the  Yp  gyro  will  cause  the 
Xp  gyro  to  dip,  thereby  ceusing  a  negative  g  component  to  be 
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sensed  by  the  Xp  accelerometer.  For  smell  t  (say  1  to  5  minutes) 
at  mid-latitudes  the  values  of  A  0X  and  A  0y  will  be  on  the 
order  of  1  minute  of'  arc.  At  these  small  values  sin  A  0X  and 
sin  A  0y  can  be  replaced  by  A  Sfx  and  A  Qfy,  respectively,  with 
negligible  loss  of  accuracy.  Therefore,  2-19  and  2-20  -can  be 
rewritten  as 


Ay  =  g*  A  0X  2-21 

Ax  =  -g*  A  tfy  2-22 

During  the  process  of  establishing  the  initial  asimuth 
wander  angle  of  a  stationary  system,  acoeleromoter  output  measure¬ 
ments  occur  only  as  a  result  of  the  system  errors  and  measurawett't 
noise.  Therefore,  the  accelerometer  outputs  Ay  and  Ax  are  a 
measure  of  the  system  errors.  Figure  2-3  is  a  block  diagram 
illustrating  equations  2-17,  18,  21,  22.  This  figure  is  a  system - 
error  model  without  the  effects  of  cross-coupling  considered. 

Cross-coupling  effects  between  the  level  axes  and  between 
the  level  and  vertical  axes  occur  due  to  platform  tilt.  The  level 
axes  cross-coupling  is  caused  by  the  level  gyros  sensing  a 
component  of  the  vertical  earth's  rate.  This  sensing  is  viewed 
from  the  navigation  frame.  For  example,  a  small  value  of  A  (2fx 
causes  the  Yp  gyro  to  lift.  This  causes  the  gyro  to  sense  a 
component  of  vertical  earth's  rate  of  magnitude  ft  zn  sin  A  0X. 

In  the  navigation  frame  the  precession  about  the  Yp  gyro  sensitive 
axis  will  be  in  the  direction  of  a  left  handed  screw.  This  cross- 
coupling  term  will  therefore  have  a  negative  sign.  For  a  small 
va-lue  of  A  the  cross-coupling  to  the  Xp  gyro  is  ftzn»sinA0y. 


Figure  2-3.  Block  Diagram  of  System  Error  Model 


without  Cross-CouDlana  Effects 


Since  A  0„  end  A  0V  will  always  be  of  small  value,  a  small  angle 
*•  / 

approximation  car.  be  used  in  lieu  of  sin  A  0  and  sin  A  0  . 

x  y 

The  level  axes  cross-coupling  is  therefore 


=  -  fl  _  -  A  0X 

zn  x 

=  ft  ,  •  A  0„ 


2-23 

2-24 


where  E  represents  the  rate  error  to  the  Y  gyro  due  to  cross- 

ycc 

co tiling  and  E*^  represents  the  gyro  cross-coupling  rate 
error.  Taking  into  account  these  level  axes  cross-coupling  terms, 
equations  2-17  and  2-18  can  be  modified  to 

A  0y(t)  =  f  f  ft  v  *  sin  0,  *  A  0„  -  ft  _  *  A  0.1  dt  2-25 


A  0y( t)  =  J  j^ft  yn*  sin  02  •  A  0Z  -  ft.^  *  A  0^jdt  2-25 
AtUt)  =  /Tn.  •  A  0..  -  ft  „  cos  0_  *  A  0_]dt  2-26 


The  remaining  significant  cross- coupling  errors  occur  as  a 

result  of  the  azimuth  gyro  sensing  the  horizontal  component  of 

earth's  rate,  as  viewed  in  the  navigation  frame.  If  the  azimuth 

wander  angle  were  zero,  that  is,  the  gyro  pointing  noi-th  in 

the  navigation  frame,  a  small  tilt  about  the  X  gyro,  A  fiL, 

P  * 

would  cause  the  azimuth  gyro  to  sense  a  component  of  level  earth's 
rate  in  the  amount  XI  yn*  sin  A  0X*  For  an  arbitrary  azimuth 
wander  angle,  the  error  rate  to  the  azimuth  gyro  due  to  level 
axes  tilt  can  be  expressed  as 

B„  =  SI  v  (cos  0Z»  sin  A  CL  -  sin  fif  •  sin  A  0V)  2-27 

-cc  Yn  •  x  z  y 

Using  the  small  angle  approximation  2-27  can  be  written 
Bgcc  «  XI  y  (cos  0Z  •  A  0X  ~  sin  0Z  *  A  0y)  2-28 

Taking  into  account  this  azimuth  gyro  error  rate,  the  azimuth 

wander  angle  error,  A  0Z,  as  a  function  of  time  is 

A0z(t)  =  J  XI  (cos  0%  •  A  0Z  -  sin  0Z'  A0y)  dt-»£0z(O) 

°  *  2-29 

Equations  2-25,  26,  29  are  used  to  construct  a  system  error 

block  diagram  which  includes  cross-coupling.  This  block  diagram 

is  illustrated  in  figure  2-4. 

To  complete  the  system  error  model  the  effects  of  inertial 
component  drifts  and  measurement  ncise  must  be  considered.  Both 
the  inertial  component  drifts  and  the  measurement  noise  can  be 
described  as  random  processes.3 

The  simplest  rcdel  for  random  gyro  drift  rate  assumes  the 
drift  to  be  a  random  constant.  For  most  applications  a  model 
such  as  this  is  not  adequate,  since,  from  empirically  obtained 
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Cross-Coupling  Effects 


gyre  drift  run3,  euto-ccrreletion  functions  of  the  exponential 
form  have  been  obtained.3 

If  Gaussian  white  noise4  is  passed  through  a  linear  first 
order  shaping  filter,  as  illustrated  In  figure  2-5,  the  output 
will  be  characterized  by  an  autocorrelation  function  of  the 
exponential  form.  In  the  figure  U  represents  zero  mean  Gaussian 
white  noise  and  €  is  the  exponentially  correlated  gyro  drife  rate 
The  i  ilationsliip  between  the  output  autocorrelation  function, 
(Oct;  ( T  ) ,  and  the  input  autocorrelation  function,  (Q  (tK  can 


t 
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Figure  2-5.  Shaping  Filter. 


be  established  using  the  following  relationship  in  the  frequency 
domain5 

=  H  (W)  H(W)<J>gu(CU)  2-30 

where  4>  ,  (a i)  and  (<*})  are  the  Fourier  transforms  of  the 

input  »t»d  output  autocorrelation  functions,.  and^g^(w) 

are  frequently  referred  to  as  the  input  and  outjnst  power  density 
spectrums.  The  term  K(W )  is  the  transform  of  the  filter  impulse 
response  or  filter  transfer  function t  and  the  bar  over  H  denotes 
the  conjugate.  The  transfer  function  for  the  Cirat  order  «b*ping 


filter  is 

H(w)  *  jS-TT/X  2~3X 

and  the  input  autocorrelation  function  for  Gaussian  white  noise  is 

8(T)  2.32 

where  (T^  is  the  variance  of  the  white  noise  end  B  (t)  is  the 
Dirac  delta  function.  The  input  power  density  spectrum  is 


UU 


(W)  = 


1*00 


*  -i»V  z 
o{t)9  dr.*cry 


2-33 


which  represents  «  constant  power  spectral  density.  Using  equation 
2-30,  the  output  power  density  spectrun  can  be  written  as 
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0  - - — r)( - • — r)0*2 

-jwtj/A  JW+  i/X  u 


or  in  the  time  domain 

-CD 

l  r 


(7aTH7X)liw+i/x  *  °fe 

Evaluating  the  above  integral®  yields 


}tur 


dw 


IT1  X  _2  -|r,/X 
(r>«  70e 


or 


2-34 


2-3S 


2-35 


2-37 


( P €€'”  2 

where  T  is  the  difference  between  two  points  in  time,  i.e.  (  t ^ 
■  t2)* 

The  v'ariance  of  the  output  random  process  is 

(0)*c r2  *  A  cr2  2-38 


4 


€€ 

Prom  empirical  data  the  standard  error,  ,  for  an  ex¬ 
ponentially  correlated  random  drift  can  be  determined.-  In 
addition,  if  an  autocorrelation  function  is  plotted,  the  time 
interval  for  the  function  to  decrease  from  <7^  2  to  a  value  of 
CT^/k  can  also  be  determined.  This  quantity,  designated  \  , 
is  aoKiatiiaea  referred  to  as  the  correlation  time.7 

If  €  y,  €x,  and  €z  are  used  to  represent  the  various 
random  gyro  drift  rates,  the  platform  tilt  angles  and  error  in 
the  aaimuth  wander  angle  due  to  these  random  drifts  (  A#€y > 
A#  €  **  and  A0£»)  can  be  computed  from 


&tf€y(s) 

* - 

11 

2-39 

J  o 

■/: 

€xdt 

2-40 

A0€s(t) 

•r 

€s  dt 

2-41 
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By  combining  these  three  equations  with  2-25,  26,  and  29  total 
expressions  for  A0V»  A0X  end  A02  cen  be  written  as 

r  * 


i: 


A0V(  t)  = 

'  cl 

&Slx(t) 

A  <$*(*)  = 


f 


ft 


yn 


ft. 


6n 


► 

1 

z 

*  A0  -  ft-  *  A0  •> 

Z  Zr  X 

€ 

y 

dt 

2-42 

: 

- 

ft  cos  0  *  A  0-  + 

”n  25  s 

€  x 

dt 

2-43 

ft  (cos  0Z*  £0  -sin  ^ 


0Z‘ 


2-44 


The  accelerometers '  outputs  due  to  random  drift  effects 
Can  be  represented ,  as  in  the  case  of  the  gyros'  random  drift 
rates,  by  linear  first  order  shaping  filters  driven  by  Gaussian 
white  noise*  These  accelerometer ‘drifts  will  be  represented  by  {X  ^ 
end  CL  y,  for  the  and  y  accelerometers  respectively.  Empirically 
derived  values  for  the  standard  error, <T^  ,  and  for  the  so-called 
correlation  time,  X  ,  can  then  be  used  in  equation  2-39  to  find 
the  variance  of  the  input  white  noise* 

Steedo®  disturbances  measured  by  the  accelerometer*  (such  ns 
vibrations,  wind  gusts,  etc.)  and  ecceleroneter  outputs  due  to 
circuit  noise  (e-g.,  amplifier  noise  in  a  force  rebalance  loop) 
can  bn  acpronimated  by  additive  white  noise  to  the  accelerometers' 
outputs.  The  actual  spectre!  density  cf  this  noise  will  be  e 
function  of  the  environment  of  the  inertial  platform,  however, 
for  convenience  an  additive  white  noise  distusbsr.ee  will  be 
assumed.  This  noise  yill  be  denoted  by  Vx  and  Vy  for  the  and 
Yp  accelerometers  respectively* 

Including  the  effects  of  the  sccaisst'Jtater's  renders  drifts 
or  id  the  additive  noise,  equations  2-21  end  2-32  cen  be  rewritten  os 
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Ay  a  g  *  A  +  ay  +  vy  2-45 

A„  =  ~g  *  A  0  +  Q  +  V  2-46 

*  y  x  x 

A  complete  systera  error  model  can  new  be  constructed  frc*. 
equations  2-42  through  2-46.  A  block  diagram  t>£  this  model  is 
Illustrated  in  figure  2-6. 
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III.  Formulation  of  State  Space  Equations 
The  system  error  model  illustrated  in  figure  2-6  is  a  multiple 
input-multiple  output  linear  system.  In  recent  years  matrix 
equations  have  been  applied  extensively  in  the  analysis  and  syn¬ 
thesis  of  multivariable  systems. 1,2  In  this  chapter  the  system 
error  model  developed  in  chapter  II  will  be  described  by  two 
matrix  equations.  These  equations  are 


X  = 

AX  + 

U 

3-1 

Y  = 

MX  + 

V 

3-2 

In  equation  3-1,  A  is  the  system  dynamics  matrix,  X  is  the  error 
state  vector  (bar  under  letter  denotes  vector),  and  U  is  the  input 
forcing  function.  A  dot  over  the  vector  indicates  differentiation 
of  the  vector  with  respect  to  time.  In  equation  3-2;  Y  is  the 
measurement  vector,  M  is  the  measurement  matrix,  and  V  is  the 
additive  measurement  noise.  Equation  3-1  is  frequently  referred 
to  as  the  system  dynamics  equation  and  equation  3-2  as  the  system 
measurement  equation. 

Equations  2-42  through  2-44  describe  the  dynamics  of  the 
system  error  model.  Differentiating  these  equations  with  respect 
to  time  yields  the  following  first  order  differential  equations 
A  iy  =  SI  •  sin  0Z  •  A  0Z  -  •  A  0X  +  €y  3-3 

Ad,  =  •  A  S(y  -  fl  ^  •  oca  ^  •  A  .  «x  3-4 

A  #z  =  (cos  0Z*  A  0X  -  sin  0Z  •  A  0y)  +  €  x  3-5 

For  the  moment,  letting  the  state  vector  X  be  composed  of 
the  elements  A  0y,  A  0X,  and  A  0z;  then  equations  3-3,  4,  and  5 
can  be  combined  to  form  the  following  first  order  matrix  differential 
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The  forcing  function  vfotor  in  equation  3-6  is  composed  of 
the  various  gyro  random  drifts.  In  chapter  II  these  drifts  were 
modeled  by  a  first  order  linear  filter  excited  by  zero  mean 
Oattssian  white  noise  (see  figure  2-5).  Therefore,  first  order 
differential  equations  describing  the  dynamical  behavior  of  €  , 
€  ,  and  €  are 


-  1/X  €  * 

y  y  y  y 

€  =  -  1/  X  €w  +  u 

x  xxx 


3-7 
3-8 

[  A  * 

€m  s  -  l/\  €  +  u  3-9 

*  as  z  z 

where  the  y,  x,  and  z  subscripts  are  used  to  identify  the  X  and  U 
associated  with  each  gyro  random  drift.  Equations  3-7,  8,  and  9 
can  be  combined  to  form  the  matrix  differential  equation 
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Similarly,  the  random  drift  of  the  accelerometers  can  be 
described  by 

°y  =  "  l/^y  ay  +  wy  3-li 

QX  =  -  l7Mx  <*x  +  vx  3-12 

where  fl  is  analogous  to  the  X  in  the  gyro  drift  model  and  W  is 
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analogous  to  U.  In  matrix  form  tha  dynamics  of  the  random  drift 
of  the  accelerometers  can  be  written  as 


The  dynamics  of  the  system  error  modal  (inclusive  of  the 

inertial  drift  models)  are  described  by  equations  3-6,  -3-10,  and 

3-13.  If  the  system  error  state  vector  is  defined  as 

A0 

y 

As!x 

AH. 


X  = 


3-14 


then  an  equation  in  the  form  of  equation  3-1  can  be  formulated  for 
the  total  system  error  model.  This  first  order  matrix  differential 
equation  is 
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Thu  accelerometer  measurements  due  to  system  errors  are 
(from  chapter  .) 

Ay  =  +  ay  +  Vy  3-16 

A,j  =  -g‘A  efy  a  x  +  Vx  3-17 

If  the  measurement  vector  is  now  defined  as 


Y  = 


y 


3-16 


tha  measurement  equation  can  be  written  in  the  form  of  equation 
3-2  as 


*x 


OgOOOOlO 
-g  0000001 


A0y 


3-19 


Equations  3-15  and  3-19  are  the  state  space  equations  for 
the  error  model  illustrated  in  figure  2-5.  Examining  the  state- 
space  dynamical  equation  (3-15)  reveals  that  the  A  matrix,  some¬ 
times  referred  to  as  the  system's  dynamics  matrix,  certains 
elements  which  ere  functions  of  0Z,  the  unknown  azimuth  wander 
angle.  Therefore,  the  angle  to  be  determined  is  e  parameter  of 
the  system.  Problems  of  this  nature  are  frequently  referred  to 
es  parameter  estimation  problems. 


''■■rVBrgXTfTSr'' 


rr^rr^rrT~’-<  [-r^wV-^r:" 
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IV.  Parameter  Estimation 

The  system  error  model  has  been  described  (in  cheater  III) 
by  the  following  matrix  equations 

X  a  A  (Ct)  X  *  U  4-1 

r  =  a  x  +  v  4-2 

»diers 

X  is  an  8x1  state  vector  of  the  systm  errors, 

U  is  an  8x1  Gaussian  white  noise  vector  that  is  the  input 
to  the  system  error  model, 

Y  is  a  2x1  vector  that  contains  the  outputs  or  amsaiarawott* 
of  the  system  error  model, 

A  (if*)  is  an  6x8  matrix  representing  the  dynamic*  of  the  *yst©* 
and  containing  elements  which  are  function*  of  •the 
azimuth  wander  angle, 

U  is  e  2x8  matrix  which  linearly  relates  X  to  Y  (uauftlly 
referred  to  aa  the  measurement  matrix),  end 

V  i«  a  2x1  Gaussian  white  noise  vector  representing  the 
additive  noise  on  the  measurement. 

A  matrix  block  diagram  of  the  system  error  model  i*  preempted 
in  figure  4-1.  Th*  wide  lines  indicate  vector-signal  flow,  and 
the  transfer  function  1/3  repreueats  8  integrator*  arranged  *o 
the  output  of  each  is  a  state  variable.  The  dynamics  uvatr&x 
A  (G^.)  indicates  how  the  outputs  of  the  integrators  «re  fed  b*ek 
to  the  inputs  of  the  integrators.  »or  example,  the  Ajj  element 
represents  the  transfer  function  between  the  output  of  the 
integrator  and  the  input  to  the  i**1  integrator. 


■A" n*  '‘■'~-~“-*F?g!tpg&!KF  r<'w-"’Z7-*-  "  t.i  «?-'> 


A»*1 
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Fitritra  4~1 ,  Matrix  Block  Diagram  cf  Syetcsa  .Srror  Model. 
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The  90 ml  of  this  investigation  ic  to  fortiylate  a  technique 
to  deteaeine  the  azimuth  waiter  aagls,  &z>  in  ainiaiat  fciiae.  In 
addition^  ££  the  gyro  drifts  {  C^)  and  the  accelerometer 

drifts  ({£  x,  Ct v)  can  be  estimated,  then  suitable  compensation  can 
be  implewtanttd  iu  the  plet-fcrn  controller.  Therefore,  astinates 
cf  the  state  of  tfets  syeton  errors,  jt,  end  an  unknown  parameter, 
fifr#  ere  heeded. 

The  method  for  eBtl^eting  the  state  of  a  system  developed 
by  Kaiaar. , 1  using  the  sainiscUm  variance  criterion,  can  be  applied 
to  linear  stochastic  ayst>«as  with  S$u#si«n<  statistics  where  the 
values  of  the  dynamical  per  sate  t«r3  ate  specified «  In  this 
investigation  one  cf  the  parameters  of  the  dynamic  taodal  is  icRown 
iwpreciflf^ly,  and  furthermore,,  precise  knowledge  of  the  parameter 
la  nocesaery  for  istabliAhiiyg  $>.H3  ift&tiel  conditions  of  the 
ioertiol  platform. 

One  jRethtd  fer  treating  this  'type  of  pxufclen  is  discussed 
by  Sorenson.*  Sorer ^0 's.  ,up£-£*ijuh  considers  the  initial  estiuated 
value  of  the  namutter  .-so  a  known  sons  tent.  A  conventional 
Kcjla.au  filter  is  theto  uynth.as£zftd.  where  owe  element  cf  the 
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=s«ti*etfcd  state  vector  is  related  to  thift  parameter  {or  is  the 
parameter),  The  parameter  value  is  then  updated  in  accordance 
with  the  state  estimate.  The  Kola«i  filter  provides  an  aatattate 
cf  the  parameter  as  well  as  the  variance  of  the  error  in  the 
estimate.  When  the  variance  reaches  an  a  triori  determined 
value  the  process  can  be  tarai  noted  and  e  value  for  tba  perameter 
established.  This  parameter  estimation  procedure  may  converge 
to  the  true  value  of  the  parameter,  however,  the  convergence  nuty 
net  he  optimal  in  x i*e.  J*or  the  problem  under  consideration 
rapid  convergence  is  of  prime  importance,  therefor**}  a  tint 
optical  procedure  must  be  utilicod* 

Uagill^»<  ha*  treated  the  problem  of  tine  csptlmel  parameter 
Ustix^tion  where  the  par  master  ie  coaatfeined  to  a  finite  mssber 
of  posslbl*  veluee.  SieBKSntel  linear  wa&ataatcris  (or.  Kplwttt 
filters)  ere  ooaetructed  for  each  of  t»tr  possible  parameter  values, 
end  the  cutputa  of  these  weighted  t*  fo\'sa  «n  optimal 

d»a|tt»4»t«;  of  the  atate  vector*  Magill  thho  sbovas  $2iat  given  the 
system  output*  the  optimal  ternneter  eciiaate  is  the  sum  o&’  tho 
elemental  perametera  each  wilted  by  a  conditional  pYobefei  lity« 
These  ^caditicnal  probabilities  ere  calculated  c?n  a  y*>cu£,elvG 
baeift  yntil  o3e  converges  to  a  wait*?  1  while  the  rent  converge 
to  z«ro.  Pigriira  4-2  illtsotratea  5tfagili’fs  tisw  «ptita»3.  ^stlsator 
an  3Ct  ccijld  ha  applied  to  this  jcroblew  of  rapid  detersfcinetioy  of 
the  esi&t-th  ^hinder  settle.  The  superscripts  os  the  refer  to  the 
veriv»via  vetuws  of  the?  !?ar ttxjter,  whi2<&  th>>  ©jper^cripts  on  X  jjrw 
diviate  twe  <vs*tbsa.ted  state  v*e~or  ss^ciated  with  ft  particular 


parameter  value. 


A,  *  , 
X‘  * 


pcfc'/yctn 
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*** 


Flg»r«4^g-_Tii>e  Optiaal  Parameter  Estimation  fjfaoillK 

Uapill's  method  is  useful  only  if  one  of  the  Sf^  (i  =  lt  2,  .. 
tk)  cbooem  for  the  el«*tntel  filler*  is  the  true  value  of  the  par* 
et*j%  ?o  practically  realize  this  condition  for  the  problem  under 
consideration,  vhere  0Z  any  essuze  0  continuum  of  values  between 
0  at»d  2tr,  a  l&rge  number  of  elemental  filters  would  be  required. 
J?br  oxampic,  <jivcn  en  initial  estimate  of  <?z  with  a  standard 
error  of  i.i5  degrees 5  it  will  be  necessary  to  construct  an 
elemental  filter  for  every  six  minutes  of  arc  in  the  range  <5  -  4.5 
to  0Z  +  4.5°,  assuming  that  the  desired  accuracy  is  of  the  order 
of  t>  uinutas  of  arc*.  Therefore,  it  can  he  seen  that  the  applica¬ 
tion  of  JI0g2.il  *s  parameter  estimator  tc  the  problem  at  hand  would 
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lead  to  computational  difficulties,  since  about  90  elemental 
filters  would  be  needed. 

For  the  problem  under  consideration  a  time  optimal  parameter 
estimator  which  allows  the  parameter  to  assume  any  of  a  continuum 
of  well  defined  values  is  needed.  An  estimator  cf  this  type  can 
be  synthesized  by  combining  the  method  of  Magill  with  Sorenson's 
approach.  As  in  Magill's  method,  this  parameter  estimator  will 
utilize  a  bank  of  elemental  linear  filters  each  initialized  with 
a  value  for  the  parameter.  However,  instead  of  fixing  the  param¬ 
eter  value  for  each  elemental  filter,  this  value  will  be  allowed 

to  converge  to  the  true  parameter  value  by  continuously  updating 
.i  a  * 

0  based  on  >*  ^0,  ),  as  in  Sorenson's  method.  Therefore,  the  new 

method  combines  the  time  optimal  feature  inherent  in  Magill 's 

approach,  with  the  continuous  value  feature  of  Sorenson's  method. 

A  limited  amount  of  filters  >.an  be  used  since  the  filters  have 

tiie  flexibility  to  change  the  parameter  value. 

The  next  step  involves  the  determination  of  a  weighting 

function  for  each  of  the  elemental  filers  so  an  overall  error 

state  vector  estimate  and  parameter  estimate  can  be  formed.  As 

will  be  shown  in  chapter  V,  the  measurement  vector,  Y,  is  differ- 

A 

enced  with  the  expected  measurement,  Y,  at  the  input  of  the 
elemental  filter.  This  expected  measurement  is  obtained  simply 
by  pre-multiplying  the  propagated  error  state  vector  by  the 
measurement  matrix.  This  vector  difference  will  be  denoted  by 
and  defined  as 


Y 


Y 


A 

Y 


4-3 


Ml 
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Referring  to  the  bank  of  elemental  filters  it  is  evident 

that  the  filters  with  parameter  values  close  to  the  true  value 

A 

w-~l  generate  expected  measurement  vectors,  Jf,  closer  in  value 
to  the  actual  measurement  vector,  Y,  than  those  whose  parameter 
value  ia  less  accurate.  Therefore,  the  outputs  of  the  elemental 
filters  where  Y  is  relatively  smaller  should  receive  the  greater 
weight.  A  simple  method  for  generating  weighting  functions  in 
accordance  with  the  above  criteria  first  requires  the  definition 

ft# 

of  the  norm  or  length  of  Y  as 

II ill  -  (i.  i)l/2  *  tv T  l)in  +•< 


Next,  an  expression  must  be  generated  which  weights  the 
elemental  filter  outputs  inversely  proportional  to  the  norm  of 
the  filter's  Y.  In  addition,  the  sum  of  all  the  weighting 
functions  generated  must  equal  1.  An  expression  which  conforms 
to  both  of  the  above  stipulations  is 

“(<>  ■  £  X*  il^3)li 

zPOli 

j=l 

1  t*’ 

where  Ci>  (0  )  is  the  weighting  function  for  the  ivn  elemental 

z 

filter  and  |J  Y  (0  ~3)  j|  the  norm  the  difference  between  the 

til 

expected  measurement  vector  as  generated  by  the  j  elemental 
filter  and  the  actual  measurement  vector. 


As  the  various  parameter  values  approach  the  true  value 
A  £  «*»  •  £ 

(V  -  - n),  the  vectors  Y  (0s  )  will  all  approach 

the  same  value.  Inspection  of  equati.  ■  4-5  shows  that  under  these 


cricumsttnces  the  weighting  functions  will  also  all  approach  the 
sauna  value  resulting  in  a  uniform  distribution  of  the  weighting 
functions. 

Finally,  a  measure  of  the  performance  of  the  parameter  estima 
tor  can  be  obtained  from  a  weighted  sum  of  the  variance  terms 
contained  in  each  of  the  elemental  filters.  The  same  weighting 
functions  generated  by  equation  4-5  will  be  used. 

Figure  4-3  illustrates  the  continuous  time  optimal  parameter 
estimator  described  in  this  chapter.  The  word  continuous  is  used 
in  the  sense  that  the  parameter  can  assume  any  of  a  oontinuum  of 


values. 
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V.  The  Elemental  Filter 


The  elemental  filter  for  the  parameter  estimation  technique 

presented  in  the  previous  chapter  can  be  synthesized  quite  readily 

using  the  methods  of  Kalman  and  Buoy.1*2  The  error  model  dynamics 

equation  (equation  4-1)  and  the  measurement:  equation  (equation 

4-2)  are  the  basis  for  synthesizing  an  elemental  filter  for  a 

specific  fifz.  The  system  error  model  is  linear,  continuous,  and 

has  a  white  noise  forcing  function.  The  measurement  noise  is 

also  white  and  is  additive.  Due  to  the  nature  of  the  problem  and 

airborne  computer  limitations  the  filter  must  necessarily  be 

discrete  performing  error  state  estimates  on  a  continuous  system. 

The  system  measurements  (that  is,  the  accelerometer  outputs)  will 

ba  observed  for  a  finite  period  (for  example,  60  seconds)  and  then 

these  measurements  smoothed  in  some  manner  and  fed  to  the  elemental 

A 

filters.  The  filters  will  then  estimate  the  state  vector,  X. 

The  various  estimates  of  the  state  vector  will  then  be  weighted 
in  accordance  with  equation  4-5  to  form  the  overall  estimate. 

The  parameter  estimate  will  then  be  formed  by  a  weighted  sum  of 
the  individual  parameter  values.  The  only  difference  between 
the  ordinary  Kalman  filter  form  and  the  elemental  filter  used  in 
this  parameter  estimator  is  feedback  of  a  component  of  the  output 
of  the  filter  to  the  elemental  filter  dynamics  matrix.  Specifi¬ 
cally,  the  filter  estimate  of  the  error  in  0Z  is  used  to  update 
the  value  of  02  in  the  system  dynamics  matrix  of  each  filter. 

The  elemental  filter  has  as  its  objective  the  generation  of 

A 

an  estimate  of  the  state  vector,  X,  svjh  that  the  mean  square  error 


42 


of  the  estimate  Is  minimized. 

A  matrix,  P,  celled  the  covariance  matrix  is  defined  as 


where  E  is  the  expectation  operator.  The  P  matrix  contains  the 

second  order  statistics  of  tho  system  errors.  That  is,  the 

diagonal  elements  are  the  variances  of  the  components  of  the 

estimation  error  and  the  off-diagonal  elements  are  the  covariances 

of  the  estimation  errors  taken  two  at  a  time.  Therefore,  the 

A 

filter  generates  an  X  which  minimizes  the  variances  on  the 
diagonal  of  P. 

The  basic  recursive  equations  for  a  diacrete  Kalman  filter 
are  adequately  treated  in  many  publications. *  The  following 
paragraphs  present  in  summary  form  the  equations  which  will  be 
utilized  to  form  the  discrete  elemental  filter  in  the  parameter 
estimator  of  this  investigation. 

Tho  error  model  dynamics  equation 

X  =  A  (02)  X  +  U  5-2 

has  a  solution  which  can  be  vritten  as 

X(T)  =  ♦  (T,0)  X(0)  +  G(T)  5-3 

where 

<t(T,0)  =  state  transition  matrix  over  the  interval  0  to 

T®>* 

G(T)  =  the  response  of  the  system  to  the  white  noise 
vector,  U,  over  the  interval  0  to  T. 

X(T)  =  state  at  time  T 
X(0)  =  state  at  time  0 


If  the  elements  of  the  White  noise  forcing  function ,  U,  are 


uncorrelatad,  a  diagonal  matrix  can  be  obtained  from  B 
If  this  matrix  is  denoted  by  the  letter  Q.  ah  approximate  expreatsion 
for  the  covariance  matrix t  P,  at  time  T  is4 

P(T)  =  «£(T,0)<P(0)*$T  (T,0)  *  |  Jor^CT.oHKT^)  *qJ  5-4 


where 


F(C)  =  ccvarianca  matrix  at  time  0 

P(T)  «  covariance  matrix  at  time  T 

^^(T,0)  =  transpose  of  the  transition  matrix 

At  this  point  it  becomes  necessary  to  distinguish  between  the 
A  A 

best  estimate  of  X  just  prior  to  fil  taring  t  X^T""},  ©nd  ju»t  after 
running  the  discrete  filter  X(T*).  T  will  be  used  fco  indicate 
values  just  prior  to  running  the  filter  and  T+  to  indicate  values 
just  after  the  filter  is  run.  This  terminology  will  also  apply  tc 
the  covariance  matrix,  P. 

At  time  0  the  filter  is  initialized  with  an  estimate  cf  the 
state  vector  X(0)  and  a  covariance  matrix  P(0).  £(C)  111x1  P(0)  ere 

then  propagated  to  the  Kalman  update  tine  T  using 

X(T’)  =  <X>  .T,0)-$(0)  5-S 

X(T“)  and  equation  5-4  for  the  computation  of  P(T~). 

A 

How  substituting  X(T~)  for  X  in  the  measurement  equation 


(4-2)  yields 

Y(T}  =  U  X(T~)  5.6 

A 

where  Y(T)  is  referred  to  as  the  expected  measurement.  The 


measurement  additive  noise  vector,  V,  does  not  appear  in  5-6  since 
the  best  estimate  of  white  Gaussian  noise  is  zero.  The  difference 
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> 


X 


between  the  expected  ftteasureraent ,  Y(T),  and  the  actual  measurement, 
Y(T},  will  be  denoted  by  Y(T).  That  is, 

Y(T)  -  Y<T)  -  Y/T)  5-7 


A  new  estimate  of  the  state  vector  can  now  be  formed  using 
the  equation 

ICT*}  *  X(T“)  +  KY  5-8 

where  K  is  the  commonly  referred  to  Kalman  gain  matrix.  This  K 

.matrix  is  ccsaputed  so  as  to  minimize  the  diagonal  or  variance 

terms  Of  the  propagated  covariance  matrix.  An  equation  which 

accomplishes  this  minimization  is7 

r  *1 

X  =  p(t")-mt  [m-p(t*)*mt  +  ap1  5-9  ' 

where  Ii  =  E 

After  the  new  state  estimate  is  formed  the  covariance  matrix  is 
updated  by* 

P(T+)  =  (I  -  KM)-P(T-)  5-10 


Equations  5-5  through  5-10  contain  the  necessary  recursive 
operations  to  construct  a  Kalman  filter  which  operates  at  discrete 
times  on  a  continuous  system.  To  further  generalize  these  equations 
replace  0  by  t  and  T  by  t  At  where  A  t  is  the  time  interval 
between  filter  runs.  The  number  of  filter  iterations  can  be 
controlled  by  checking  the  variance  terms  cf  tha  P  matrix. 

When  a  desired  variance  for  an  error  state  is  reached  the  problem 
con  bo  terminated. 

For  the  parameter  estimation  technique  presented  in  this 

investigation  it  is  necessary  that  the  value  of  0  (the  azimuth 

z 

wander  angle)  in  Che  A  matrix  be  updated  by  the  estimate  of  A  0_ 


cxt'fftr  each  filCer  run.  This  update  is  #cco*rpi£aU«£d  by  suaming 
the  e*t£*u*t:<n  of  t he  etfror  in  s£z  (that  is  &  <$„  (t  «■  &  t*})  with 
the  osirrsut  value  of  cJRr  Ther^fsre,  if  0Z  (t  *  &  t+)  £«  used 
to  denote  the  updated  value  of  02  and  0J,  {t  +  &  t“)  S?»e  value 
b*£asxs  updafc issg 

1  (t  +  At+)  s  0‘  (*  *  A  t“)  «  A  0  (t  ♦  A  t+)  5*12 
This  updated  valuta  of  the  azimuth  wander  angle  will  then  b©  used 
in  th*  forwation  of  a  new  A  matrix  for  the  eleaecntel  filter* 

A  block  diagram  of  the  elemental  filter  £©  presented  £n 
figure  5-1  and  an  equation  flow  chatt  for  the  filter  algorithm 
presented  in  figure  5-2. 


Pigure  5-1.  Block  Diagram  of  Eleraentel  Filter. 
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PROPAGATE  COVARIANCE  MATRIX  TO  *+Afc“  m  j 

P(  t+At“  )=<K  ti-Ac,  *  )  ♦  P(  t  )•  $  T(  t+At ,  t  )+|&jQ$f  ( t+At  *  t  )*$(  t^Ac ,  t  )oj| 

PROPAGATE  ERROR  STATE  ESTIMATE  TO  t+A*~  j 

X ( t+At“  )  =  $  (t+At,t)*$(t) 


COMPUTE  KALMAN  GAIN,  K 
K  =  p(t+At-)*uT  [up(t*At‘)  MY  +  r] 


-I 
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FORM  EXPECTED  MEASUREMENT 

Y(t+At)  =  M  x  (t+At“) 

DIFFERENCE  SYSTEM  MEASUREMENT  WITH 
EXPECTED  MEASUREMENT 

Y(t+At)  -  Y(**&t)  -  Y( t-i-At ) 

FORM  UPDATED  ERROR  STATS  ESTIMATE 
“  X(t+At')  +  K  Y(t+At) 

UPDATE  COVARIANCE  MATRIX 
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•  *  X 
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ELEMENTAL 
FILTER  ERROR 
STATE  OUTPUT 


Figure  5-2.  Elemental  Filter 
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VI.  The  System 

The  prime  function  of  the  platform  controller  during  the 
initialization  process  is  to  apply  earth  rata  torqueing  signals 
to  the  platform  gyros.  The  level  axes  torquing  rates  era  com- 
puted  using  equation  2-5  with  #2,  the  azimuth  wander  angle) 
replaced  by  its  estimate ,  $8.  The  desired  equilibrium  condition 
is  obtained  when  the  correct  azimuth  wander  angle  is  determined 
(within  a  specified  accuracy)  and  the  gyro  and  accelerometer 
random  drifts  compensated. 

Using  an  external  reference  or  from  observation  of  the 
platform  without  earth  rate  to r quin g,  the  initial  azimuth  wander 
angle  can  be  estimated  with  an  accuracy  of  1.5  degrees  (lCT).  A 
number  of  values  are  then  selected  in  the  vicinity  of  this  initial 
estimate  and  an  initial  weighting  distribution  for  these  values 
assumed.  The  most  intuitive  approach  to  selection  of  these 
initial  values  is  to  assume  a  symmetrical  distribution  with  sym¬ 
metrical  weighting.  For  example ,  if  02  is  the  initial  estimate) 

i  .  „ 

values  of  0Z  for  the  elemental  filters  could  be  fifa,  0a  ±  1°, 

0a  t  2°,  etc.,  end  an  initial  triangular  weighting  centered  on 

0mi  The  initial  or  a  priori  weighting  distribution  is  utilised 

only  for  the  first  run  of  the  parameter  estimator.  Subsequent 

runs  use  the  weighting  factors  given  by  equation  4-5, 

_  a 

The  error  state  vector,  X,  for  each  of  the  elemental  filters 
is  initialized  with  all  elements  equal  to  zero.  The  covariance 
matrix  for  each  of  the  elemental  filters  is  initialized  with 
estimated  variances  for  each  error  along  the  diagonal  and  all 


off-diagonal  terras  zero.  After  the  first  run  of  the  parameter 
estimator  each  of  the  elemental  filters  contains  an  estimate  of 


the  parameter,  0Z,  estimates  of  the  angle  misalignments  (  A  0X, 

A  0y,  A  02),  estimates  of  the  inertial  component  drifts  (f  x» 

€  ,  <  ,  fl  .  CL  ),  and  an  updated  covariance  matrix,  A  new  elemental 
y  *  *  y 

filter  parameter  estimate  is  formed  using  equation  5-12.  An 
overall  parameter  estimate  can  then  be  formed  using  the  a  priori 


weighting  factors  in  accordance  with 

n 

■X  ■  Z 
i=l 


a 


6-1 


where  n  ia  the  number  of  elemental  filters.  This  value  of  0Z  is 
then  used  in  the  computation  of  the  earth  rate  torquing  signals 
generated  by  the  controller. 

Compensation  for  the  inertial  components'  drifts  is  accom¬ 
plished  by  weighting  and  summing  each  of  the  elemental  filter 
drift  estimates,  and  then  applying  the  results  directly  to  the 
inertial  components.  The  inertial  components'  drift  corrections 
are  then  subtracted  from  the  inertial  component  drift  estimates 
in  each  of  the  elemental  filters.  This  leaves  a  residual  value 
in  the  fourth  through  eighth  positions  of  the  error  state  vector. 
This  residual  value  will  be  propagated  to  the  next  update  time 
in  accordance  with  the  propagation  equation  contained  in  the 
elemental  filter  algorithm.  The  A  0X  pnd  A  &y  angle  misalignments 
are  torqued  to  zero  by  the  leveling  loops  and  the  corresponding 
error  states  in  each  of  the  elemental  filters  set  to  zero.  The 


A  02  term  in  each  of  the  elemental  filters  is  set  equal  to  the 
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difference  between  the  overall  estimate  of  0Z  and  the  value  of 
0B  in  the  particular  elemental  filter. 

On  the  second  and  subsequent  parameter  estimator  runs 
weighting  factors  are  computed  in  accordance  with  equation  4-5. 
These  weighting  factors  are  then  used  in  the  computation  of  the 
overall  parameter  estimate  and  the  inertial  components *  compensa¬ 
tion. 

The  primary  result  is  to  obtain  a  value  for  the  azimuth 

wander  angle. within  a  specified  statistical  variance.  A  variance 

for  the  parameter  estimator  is  formed  by  multiplying  the-  weighting 

factor  for  each  elemental  filter  by  the  variance  of  A  0„  in 

the  elemental  filter  covariance  matrix.  When  this  variance  falls 

balow  a  specified  value  the  initialization  procedure  is  terminated. 

The  current  value  of  0  is  maintained  in  the  controller  earth 

z 

rate  computation  and  the  current  inertial  components'  compensation 
maintained.  If  in  the  navigation  an  augmented  system  with 

a  Kalman  filter  is  utilized,  the  azimuth  wander  angle  and  the 
inertial  components'  compensation  will  continue  to  be  corrected 
during  flight. 

Figure  6-1  is  a  block  diagram  of  the  rapid  initialization 
technique  described  in  this  investigation. 


j£2  Bl.PPk  Diegraa  of  Rapid  Initialisation  Tgchni 
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VII.  Rapid  Initialization  Simulation 

‘/To  verify  that  the  parameter  estimation,  technique  described 
in  this  investigation  vill  converge  to  the  true  value  of  the 
parameter  (to  a  statistical  accuracy),  a  computer  simulation  of 
the  entire  system  was  programmed  for  use  on  a  Burroughs  -  5500 
computer  using  Extended  ALGOL.  The  inertial  platform  was  simu¬ 
lated  by  first  order  difference  equations  which  were  iterated 
every  second  to  approach  a  continuous  system.  The  outputs  of  the 
platform  accelerometers  were  smoothed  using  straight  averaging 
and  fed  to  the  filters.  The  recursive  elemental  filter  equations 
were  contained  in  an  ALGOL  Procedure  named  Kalman  which,  when 
called  up,  ran  through  ths  equations  delineated  in  figure  5-2 
using  the  appropriate  $,  and  P  for  the  particular  elemental 
filter.  The  norm  of  Y  was  also  computed  in  each  filtar  to 
facilitate  computation  of  the  weighting  factors.  Ths  filter  was 
run  every  60  seoonds.  The  covariance  propagation  equation  and 
state  vector  propagation  equation  were  run  every  10  seconds. 

This  was  necessary  due  to  the  error  introduced  by  truncation  of 
the  series  approximation  for  the  transition  matrix  at  the  second 
term..  During  the  time  the  filter  was  running,  the  platform- wss 
releveled  (approximately  5  seconds  real  time  to  relevel). 

To  oompare  this  rapid  initialization  technique  with  present 
rtate-cf-the-art  gyrocompassing,  a  thermally  compensated  experi¬ 
mental  platform  was  gyrocompassed  from  an  azimuth  error  of  1.5° 
and  the  dynamic  azimuth  error  recorded  on  a  strip  chart.  Pigure 
7-1  illustrates  the  amount  of  time  necessary  to  achieve  an 
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alignment  accuracy  of  approximately  6  minutes  of  ere  for  this 
system  (8  minutes  real  time).  Statistical  data  on  the  performance 
of  the  inertial  components  was  then  obtained  from  the  manufacturer 
and  tasad  in  the  rapid  initialization  computer  simulation.  The 
data  obtained  was  as  follows: 

Gyro  Random  Drift  (Correlated):  1(T=  .0 2°/fcR 

Correlation  time  »  6,900  men 
AocsIcvomo  tar  Random  Drift  (Correlated):  1<T~  5  x  10“  *g 

Correlation  time  *  30,000 

A  measurement  noise  vector  was  assisted  which  consisted  of  am  e>- 
oorrel'ated  sequence  which  changed  value  every  second.  Computer 
rune  were  then  made  for  azimuth  wander  angles  of  0,  30,  and  49 
degress.  Five  runs  were  made  for  smah  ?£  these  aslmnth  sngls«j 
•s$t  with  a  different  random  sequence.  The  various  rondos; 

—games  »  were  obtained  by  changing  the  6  digit  key  in  this  Algol 
Procedure  "  Independent  Gaussian  Random  Variables"  (JGKV) .  The 
parameter  eatimator  contained  5  slementol  filters  which  were 
initialised  at  values  of  *1  end  ±2  degrees  about  an  initial 
estimate  which  was  in  error  by  1.5  degrees.  An  initial  ttteqnldt 
weighting  distribution  centered  on  -the  initial  estimata  wes  ' 
assumed . 

•  The  computed  variance  of  the  parameter  estimator  indicated 
that  an  accuracy  of  6  minutes  of  arc  (1<T*)’  could  be  reached  in 
3  minutes  of  running  time  (plus  approximately  3  times  5  seconds 
for  the  filter  runs).  The  actual  average  error  of  ths  fifteen 


m 
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rurvn  after  3  minutes  of  running  time  was  5,S  minutes  of  arc. 

Table  7-1  lists  the  final  azimuth  wander  angle  errors  for  each 
of  the  runs.  Figure' 7-2  is  a  plot  of  the  average  error  for  the 
fifteen  runs  versus  time.  Typical  computed  error  variances  of 
the  parameter  estimator  as  a  function  of  time  are  illustrated  in 
figure  7-3,  and  figure  7-4  shows  the  initial  weighting  distribu¬ 
tion  and  «  typical  distribution  after  the  third  iteration. 

The  simulation  program,  is  contained  in  appendix  A  mod  the 
printout  for  three  of  the  fifteen  runs  is  contained  in  appendix  B. 


A«&autii  Srrrsr  Minutas  of  Aro 


TABLE  7-1 

Pinal  Azinuth  Wander  Angle  Error  in  Minutes  of  Arc 
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VIII .  Conclusions 

The  problem  of  rapid  determination  of  the  angular  orienta¬ 
tion  of  an  inertial  platform  with  respect  to  the  rotating  earth 
reference  frame  has  been  investigated  and  a  new  method  for 
initialization  presented.  A  linear  error  model  for  an  inertial 
platform  in  an  azimuth  wander  mechanization  woe  derived  and  the 
equations  expressed  in  state  space  form.  The  problem  in 

initialization  wee  identified  aa  a  parameter  estimation  problem. 

The  two  methods  of  solving  parameter  estimation  problems 
currently  available  in  the  literature  were  then  presented.  The 
first  method  allowed  parameter  estimation  for  a  parameter  which 
oan  assume  a  continuum  of  values;  however,  the  method  was  not 
time  optimal.  The  seoond  method  examined  was  time  optimal;  how¬ 
ever,  it  was  constrained  to  problems  where  the  parameter  is 
limited  to  e  finite  number  of  possible  values. 

By  combining  features  of  both  methods  and  using  a  simple 
weighting  scheme,  the  time  optimal  method  was  extended  to  permit 
parameter  estimation  for  a  parameter  that  can  ansosm  m  oontinmm 
of  values*  In  addition,  the  new  parameter  estimator  provides  a 
variance  term  for  the  parameter  estimate «  This  allows  the 
initialization  procedure  to  be  terminated  when  a  predetermined 
variance  is  achieved. 

The  necessary  algorithms  for  the  parameter  estimator  were 
develop  id  and  a  computer  simulation  of  the  system  performed. 
Results  of  the  simulation  show  a  twofold  improvement  in 
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initialisation  time  over  a  state-of-the-art  experimental  system 
presently  tinder  evaluation. 
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4PFSNDXX  A 

Computer  Simulation  Program 
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THIS  PROGRAM  IS  BAPl 

A  SIMULATION  OF  A*  RAPID  INITIALIZATION  TECHNIQUE 


TRUE  AZIMUTH  AIvGLECPHIZ)  *  0.0.00 

AN  ESTIMATE  OF  PHIZ  IS  MADE  EVERY  60  SECONDS- 

THE  INITIAL  RANDOM ‘ NUMBER  GENERATOR  KEY  ^  23232323 


TIME  IN  MINUTES  **  ,0.00  n 

;  •"  a.  • . ••  '  • 

THE  ELEMENTAL  FILTER  PHIZ  ESTIMATES  ARE  ' 
PHlZHAtl  *  -O.'SOO"'  »  ”  r  "  ‘  *  “ 

PHIZKAT2  £  "OV500 

PHI2HAT3  i  1Y500-  . 

PHIZHAT4  £  2Y500  '  ' ' 

PHIZHAT5  £  3V500*  .  .  : 


-THE  WEIGHTING  FACTORS  ARE  •* 
Vl"=  0.1000  r  ' 

«2  £  0Y20C0 
¥3  £  0Y4000 
«4  £  0V2000 
¥5  £  0 VI 000 


Tins  DASARO  PARAMETER  ESTIMATE^  IS  1  .5000 
THE  ERROR  VARIANCE  IS  *  $«80000ft-&4 


THE  ERROR  IN  MINUTES  OF  ARC  IS  -90.00000 

V 

TIME  IN  MINUTES  *  1.00  1 


THE  ELEMENTAL  FILTER  PHIZ  ESTIMATES ’ ARE 

PKIZIiATl  =  -0.'696 .  “  <**’  “ 

PHI2HAT2  =  -0-09S 

PHIZHAT3  £  "O'** 50  7  - 

PHIZHATA  £  I'.’IOS  . 

PHIZHAT5  =  IV?  11 
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THE  WEIGHTING  FACTORS  ARE  ' 

ui~*-a;iooo  -  -  •  *  . 

VP.  5  0V2000 
V3  *  0V4000 
s4  *  ovsooo 

k'5  *  O  VIOOO  _  . 


•  ,  » 

THE  SASARO  PARAMETER  ESTIMATE  IS  0-5072 
'  THE  ERROR  VARIANCE  IS  .2*706810-04 

i  ■ 

•  ■  .  e  *  # 

■  .  .  .  V 

THE  ERROR  IN  MINUTES  OF  ARC  IS  -30.43337 


-TIME  IN  MINUTES  =  2.00 

-  -  .  '  ‘  v 


TKE  ELEMENTAL  FILTER  PHIZ  ESTIMATES  ARE 
PHIZHAT1"*  '  -0 .'480"  •"  "  '  '  ~ 

PHIZHAT2  i  -0V071  ‘ 


...?  HI2HAT3  = 
PUX3HAT4  '= 


'0V338 

0V748 


PKI2HAT5  *  1V1S8 

mm  —  «  ».  «  rf*  mm  *»  mm. 


.  THE  WEIGHTING  FACTORS  ARE 

«r»  0 •  1 769 . •  •  •  - 

WS  *’  0V21‘80 
U3  ®  0T2391 
H4  *  0V2016 
H5  *  OVi  623 


THE  DAS  AH  0  PARAMETER  ESTIMATE  IS  0*3183 
THE  ERROR  VARIANCE  IS.  3*148060-06 

i- 

THE  ERROR  IN  MINUTES  OF  ARC  IS  '  -19*09694 
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:  .  ,  •/  -  ' 

TIME  IN  MINUTES  »  3.00  *  '  /-'**;• 

v  ’  f‘  ‘  -  •  .  •  •  .. 

t  *  j  »  ‘  • 

THE  ELEMENTAL’  FILTER  PHIZ  ESTIMATES  ARE- 
PHIZHAT1’*  * -0 •  652 7  *  *  *  "  '  " 

PHIZHATS.  *  .  HV.3'66 
•  PH12HAT3  *  '-'0V030  * 

,  FB1ZHA74  *  ‘  0V208 

.PH12KAT5  *  07497  ’  '  ;  •  * 


ai 

THIS  PRGGRAH  IS'KAPl 

A  SIMULATION  OF  A  RAPID  INITIALIZATION  TECHNIQUE 

•  i  »  .  , 

TRUE  AZIMUTH  ANGLECPKIZ)  =  30 <.000  *  ■ 

AN  ESTIMATE  OF  PHIZ  IS  MADE  EVERY  60  SECONDS 
THE  INITIAL  RANDOM  NUMBER  GENERATOR  KEY  •-  23238323 


TINE  IN  MINUTES  =  0.00 

•  •  '  •• 


THE  ELEMENTAL  FILTER  PHIZ  ESTIMATES  ARE 

PHIZHATI 29*500"  . 

PK1ZHAT2  *  30V500 

PHIZHAT3  =  31Y500 

PHIZHA'M  =  32V  500 

PHIZHAT5  =  33VS00 


1HE  WEIGHTING  FACTORS  ARE 

Wl"*  0 ilOOO . r 

V2  S-0V20D0 
W3  =  0V4000 
V4  =  0V8000 
W5  =  071000 


THE  DASARO  PARAMETER  ESTIMATE  IS  31.5000 

»  mm  rn.  ~  • 

THE  ERROR  VARIANCE  IS  6.600000-0*1 


THE  ERROR  IN  MINUTES  OF  ARC  IS  -90.00000 

TIME  IN  MINUTES  =  l  .00 


THE  ELEMENTAL  FILTER  PHIZ 
PHIZKAT1  =  29. '297 

PHIZHA72  =  29**699 

PHIZHAT3  =  30-' 501 


PHlZHATA  =  3IV10A 
FHIZHATS  -  3 1 V 70  7 


Estimates  are 


Q2 


*<*>«'  '%*•  »  **  *  .}<•*  •i's  / ;*' * ■  ,  *  ’  *  ,*  •  '■■/  *  •»*'* 

THE  WEIGHTING  FACTORS  ARE 
VI' «=  CJilOOO  ' 

W2  =  OY2O0O 
%  W3  =■  0V4000 

W4  =  GV2000 
’  W5  =  0 VI 000 


k  »*«**-  *  •  - 


THE  DASARO  PARAMETER  ESTIMATE  IS  30 .$014 
THE  ERROR  VARIANCE  IS  *  2.706810-04 


•'THE  ERROR  .IN  MINUTES  OF  ARC  IS  -30.08420 


;  .  <•’  .  ,  v  • 

«  ’  .  :  1  -■  .  *  . ...  . 

TIME  IN  MINUTES  *='  2 •00 


'  '  \  / 


THE  ELEMENTAL  FILTER  PHIZ  ESTIMATES  ARE 
.  PHlZRAtl ' 29*620“  '  *  *** 

PH1ZHAT2  &  30V030 

1 PHIZHAT3  »  30Y442 

7  PHI2HAT-4  «  30V853 
.  PH1ZHAT5  =  31Y265. 


THE  WEIGHTING  FACTORS  ARE 

wr=  q;i?37 . *  '  - 

U2  =  0V2122 

W3  =  0V2393  •  ^ 

W4  =  0*206? 

W5  =  0V1681  i  ^ 


THE  DASARO  PARAMETER  ESTIMATE  IS  .30.4351 

»  \ 

THE  ERROR  VARIANCE  IS  3^145060-06 


THE  ERROR  IN  MINUTES  OF  ARC  IS  -26.10787 


"JMWW 


.  . —  •  '1- 


■t 
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TIME  IN  MINUTES  =  3.00 


THE  ELEMENTAL  FILTER  PHIZ  ESTIMATES  ARE 

PHIZftATl"*  '29V468'  7 . 

PHIZHAT2  S  29V 755 
PHIZHAT3  «  30  VO  42 

PHIZHAT4  s  30V330 
PHIZHAT5  =  30  V  620 


THE  WEIGHTING  FACTORS  ARE 
.  Wl 0.2046  '  “  *  * 

¥2  «  0V2379 
W3  «  0V2213 
W4  =  0  Vi  860 
W5  £  0  VI 502 


THE  DASARO  PARAMETER  ESTIMATE  IS  29.9967 
THE  ERROR  VARIANCE  IS  2.983178-06 

...  w  ••  ,•  •  •  »  «  *•  •  ..  *f  *  •  •  •  «  *• 


V 

THE  ERROR  IN  MINUTES  OF  ARC  IS  0.19943 

/ 

VARIANCE  BELOW  VALUE  CORRESPONDING  TO 
.A*  l  SXGM  ERROR  OF  6'  MINUTES  '  OF'ARC 

INITIALIZATION  COMPLETE 
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THIS  PROGRAM  'iS  RAP1  \ 

A  SIMULATION  OF  A  RAPID  INITIALIZATION  TECHNIQUE 


TRUE  AZIMUTH  ANGLECPHIZ  >  =  45.000 

AN  ESTIMATE.  OF  PHIZ  IS  MADE  EVERY  60  SECONDS 

THE  INITIAL  RANDOM  NUMBER  GENERATOR  KEY  =  23232323 


TIME  IN  MINUTES  =  0.00 


THE  ELEMENTAL  FILTER  PHIZ  ESTIMATES  ARE 
PHIZKAT1  =  ’44*500' 

PHIZHAT2  =  45*500 

PHIZHAT3  =  467500 

PHIZHAT4  =  47V500. 

PH2ZHAT5  =  48  V 500  : 


i 


THE  WEIGHTING  FACTORS  ARE 

wi  "*»  o.iooo  '  "  ;  ’  ” 

W2  =  078000 
Vi 3  =-  074000 
W4  =  078000 
Vi  5  =  0  VI 000 


THE  DASARO  PARAMETER  ESTIMATE  IS  46*5000 
THE  ERROR  VARIANCE  IS  6i80000@-04  „■ 

f  •  .  *  •/ 


THE  ERROR  IN  MINUTES  OF  ARC  IS  ”90*00000 


TIME  IN  MINUTES  =  i .00 


THE  ELEMENTAL  FILTER  PHIZ  ESTIMATES-  ABE 
PKIZKAT1  =  44  •  30 1  " 

r  PHIZHAT2  -  447904 

PHIZHAT3  =  45'*  50  7 

PHIZHAT4  *  46"*"  1 10 

PHZZHAT5  =  467713 


THE  WEIGHTING  FACTORS  ARE 

VI' «  o.iooo  . 

VS  <?  OVS0OO 
W3  *  074000 
V4  =  0V2000 
VS  «  OViOOO 

•  -  .  ...  .  4 


THE  DASAHO  PARAMETER  ESTIMATE  IS  45*5068 
THE  ERROR  VARIANCE  IS  2. 70681 §-04 


t 

THE  ERROR  IN  MINUTES  OF  ARC  IS  -30.40692 


TIME  IN  MINUTES  =  2.00 

•  * 

t  *  + 

THE  ELEMENTAL'  FILTER  PHIZ  .ESTIMATES  ARE 
PHIZHAT1'=  44.688"  '  "  * 

PHIZHAT2  «  4SV099  • 

PH1ZHAT3  *  4SV5H 
PHIZHAT4  »  45V923 

PHIZHAT5  =  46V336 


THE'  WEIGHTING  FACTORS  ARE 
VI  '*  0.1706  *  "  ' 

V2  *  OV2090 
V3  = 1 072392 
W4  »  0V2098 
VS  -4  0V1714 

THE  DASARO  PARAMETER  ESTIMATE  IS  45.5123 

» 

THE  ERROR  VARIANCE  IS  3.145060-06 

...  '  . .  .  P  ... 


THE  ERROR, IN  MINUTES  OF  ARC  IS  -30.74099 
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TIME  IK  MINUTES  =  3.00 


THE  ELEMENTAL  FILTER  PHlZ  ESTIMATES  ARE  • 

PHIZHAT1  "*®  '44i 557*  \ . 

PHIZHAT2  =  44V844  ,  ' 

PHIZHAT3  *■  -4SV13I  '*  ■- 

Ptil ZHAT4  «  45V419 

PHIZHAT5  »  35V 709  v  '  '  '• 


t* 


v.  . 

'  \  ■■ ‘  ’  ••  %  ‘  ! 

the  weighting  factors  are  •*  '• 

wr*  o;so27  '  '  '  *  '■  .  •: 

¥2  ■  0V&349  •-*  >r'"\ 

¥3  *  0V2222  .  , ,  •  •-  V  -v  :.v  •  ’ 

•W4  *  0 VI 878  '  -  ■  •  ,*>;V  <  h  \  "•  ; 

¥5  «  OVi  524  :  *•':  '  •  •'  ■*'*  'l  v.' 

r-*  \s  ••  . 

.  *  ■  *  i-  .  .  •  •  •  •  . 


•  ■  -  •  ■  * .  '  i* ■  •  -  . v 

,  -  %  K*  '  •  •  'll.  * 

*  v..  .  A  a  ‘  ■  j 

■  .«  .  >  . .  .  •  **'■■•.■*  , ;  %  •  ’  r  ,  %• 


THE'  EASARO  PARAMETER  ESTIMATE  IS  45  .6894  "  J 


THE  ERROR  VARIANCE'  IS  2*9831 69-06'.  N 


THE  ERROR  IN  MINUTES  OF  ARC  IS  -5*38467 


VARIANCE  BELOW  VALUE  CORRESPONDING  TO 
N  A  1  SIGMA  ERROR  OF  6 * MINUTES 'OF'ARC 

r>  ‘  *•  •  ' 


INITIALIZATI ON . COMPLETE 
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Glossary 

local  level  right  handed  tried 0  X^  is  level 
and  pointing  east,  Y^  is  level  end  north,  and 
2n  is  up. 

mutually  orthogonal  platform  axes, 
earth's  rotation  rate  with  respect  to  the 
inertial  or  fixed  frame. 

level  component  of  earth's  rate  along  Ya. 
level  component  of  earth's  rate  along  Xn. 
vertical  component  of  earth's  rate  along  Zn. 
latitude. 

azimuth  wander  angle, 
estimate  of  <&z. 
difference  between  0'z  and 

eagles  between  the  locally  level  plane  and  the 

platform's  Y  and  X_  axes,  i.e.,  A  0„  is  the 
p  P  x 

angle  between  Yp  and  the  locally  level  plane 
measured  by  a  right  handed  rotation  about  Xp, 
similarly  for  A 

errors  in  earth's  rate  torquing  to  the  Y  and 
X  gyros  respectively  due  to  error  in  azimuth 
wander  angle. 

errors  in  torquing  to  the  Y,  X,  Z  gyros  due 
to  cross-coupling  effects. 

the  gravity  vector;  defined  to  be  positive 
in  the  vertical  upward  direction. 


g 


!x> 


Ay,  Ax 

u 


platform  level  axes  accelerometer  measurements. 


^€(r) 

,!{AU(T) 
8  (T) 

♦uo(w) 

a* 


fix*  fiy 


A 


U 

Y 
M 

Y 

A 

Y 


zero  mean  Gaussian  white  noise, 
exponentially  correlated  rondom  gyro  drifts, 
autocorrelation  function  of  random  gyro  drift, 
autocorrelation  function  of  Gaussian  white  noise. 
Dirac  delta  function. 

power  density  spectrum  of  random  gyro  drift, 
power  density  spectrum  of  Gaussian  white  noise, 
standard  deviation  of  a  random  process, 
variance  of  a  random  process. 

"correlation  times"  for  the  random  gyro  drifts} 
values  of  for  which  the  exponential  autocorrela¬ 
tion  functions  decrease  to  l/e(Ta. 
natural  logorithm  base. 

"correlation  times"  for  the  random  accelerometer 
drifts, 

accelerometer's  random  drifts. 

error  state  vector. 

system  dynamics  matrix. 

zero  mean  Gaussian  white  noise  vector. 

measurement  vector. 

measurement  matrix. 

additive  measurement  noise. 

,  i 

elemental  filter  with  02  parameter  estimate, 
estimate  of  X. 


estimate  of  Y 


difference  between  actual  measurement  Y,  end  astirae** 
A 

ted  measurement,  Y. 

#V 

norm  of  Y. 

weighting  factor  for  outpucs  of  ^iemennal  filter 
with  0„1  parameter  estimate, 

covariance  matrix  of  errors  in  state,  vector  estimate 
expectation  operator. 

state  transition  matrix  over  the  interval  t  to  t  + 

At. 

matrix  of  variances  of  input  forcing  function,  U. 

Kalman  gain  matrix. 

matrix  of  additive  noise  vstianceo. 


Index 


Accuracy  of  self-contained  alignment 

Additive  noise ,  22 

Approximations, 

moll  angle,  14,  17 
Auto -'Correlation  functions,  19 
Asicuth,  alignment,  3 
gyre,  3 
references,  3 
wander  angle e  5,  11 
Controller, 
platform,  46 
Correlation  time,  21 
Covariance  matrix,  42 
propagated ,  42 
Cross-coupling  errors,  16 
Dirac  delta  function,  20 
Earth's  rate,  10 
Elemental  filter,  41 
block  diagram,  45 
computer  algorithm,  46 
Error  model, 

block  diagram,  24 
Error  vector,  28 
Feedback, 

in  elemental  filter,  42 
Gaussian  white  noise,  19 
Gyro  compassing ,  3 
in  polar  regions,  5 
Gravity,  15 

Inertial  Component,  warm-up 
behavior,  6 
random  driftc,  6,  18 
Inertial  franc,  2 
Inertial  measurement  unit,  1 
Inertial  navigation  system,  1 
Inertial  platfors,  2 
Initial  conditions,  1 
Kalman,  R.E.,  33 
gain  matrix,  4.4 
Latitude,  10 
Leveling,  platform,  3 
Local  level  triad,  2 
Mayill,  D.T.,  34 
Magnetic,  compass,  6,  14 
variation,  14 
Measurement,  vector,  28 
matrix,  32 
expected ,  43 

Navigation  reference  frame,  2 
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Nciea,  jssamireaent,  22,  5 
Noise  vector,  32 
Norm,  of  a  vector,  37 
Paraae ter  estimation,  29,  32 
Bower  density  specttum,  20 
Rand  cat,  ecceleroneter  drifts,  22 
gyro  drifts,  31 
Random  processes,  13 
Rotation  matrix,  12 
Schuler- tuned,  2 
Shaping  filter,  20 
Simulation, 
computer,  52 
Sorenson,  H.V.9  33 
Standard  error,  21 
State  space  equations,  26 
State  vector,  23 
Stationary  system,  3 
Theodolites,  electronic,  3 
Thermal  compensation,  6 
Time  optimal  parameter  estimator,  36 
Torquing,  arrcr,  13 
Transition  Matrix,  42 
Variance, 

of  white  noise,  20 
of  parameter  astimate  drror,  33 
ttritfhting  factor,  37 
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it.  Xostnact  - 

~  The  accuracy  of  an  aircraft  inertial  navigation  system  depends  upon  the  accuracy 
with  which  the  system  is  initially  aligned.  One  procedure  for  initial  alignment  in¬ 
volves  the  use  of  an  external  reference.  This  method  utilizes  equipment  which  is  much 
too  elaborate  for  normal  operational  use.  An  alternate  procedure  uses  the  system's 
inertial  sensors  in  a  self -contained  method.  If  sufficient  time  were  available,  the 
self-contained  method  could  achieve  accuracies  commensurate  with  the  sensor  accuracies; 
however,  in  an  operational  environment  it  is  usually  necessary  to  sacrifice  some 
accuracy  in  the  interest  of  achieving  a  more  rapid  initiation.  This  dissertation  in¬ 
vestigates  the  methods  presently  available  fer  initialization  of  an  inertial  platform 
in  an  azimuth  wander  or  free  azimuth  instrumentation  and  presents  a  new  method  for 
rapid  initialization.  The  paramount  problem  is  the  determination  of  the  initial  azimuth 
angle  in  minimum  time  in  the  presence  of  random  gyro  drifts,  random  accelerometer 
drifts,  and  measurement  noise. 


